Library import
library(SpatialExperiment)
library(data.table)
library(scater) # it imports scuttle
library(cowplot)
library(ggplot2)
theme_set(theme_bw())
library(matrixStats)
library(dplyr)
library(tidyr)
library(tibble)
library(arrow)
library(scales)
library(sf)
library(tmap)
Read-in function
readCosmxSPE <- function(dirname = dirname,
countmatfpattern = "exprMat_file.csv",
metadatafpattern = "metadata_file.csv",
coord_names = c("CenterX_global_px",
"CenterY_global_px")){
countmat_file <- file.path(dirname, list.files(dirname, countmatfpattern))
metadata_file <- file.path(dirname, list.files(dirname, metadatafpattern))
# Read in
countmat <- data.table::fread(countmat_file)
metadata <- data.table::fread(metadata_file)
# Count matrix
counts <- merge(countmat, metadata[, c("fov", "cell_ID")])
counts <- subset(counts, select = -c(fov, cell_ID))
cell_codes <- rownames(counts)
features <- colnames(counts)
counts <- t(counts)
rownames(counts) <- features
colnames(counts) <- cell_codes
# rowData (does not exist)
# colData
colData <- merge(metadata, countmat[, c("fov", "cell_ID")])
spe <- SpatialExperiment::SpatialExperiment(
assays = list(counts = counts),
# rowData = rowData,
colData = colData,
spatialCoordsNames = coord_names
)
return(spe)
}
Setting directories
dirname <- "/NAS06/work/Spatial_omics/DBKero/CosMx_Breast/CosMX_data_Case2"
count_pattern <- "Run5810_Case2_exprMat_file.csv"
meta_pattern <- "Run5810_Case2_metadata_file.csv"
sample_name <- "CosMx_DBKero_BC"
SpatialExperiment object creation and exploration
spe <- readCosmxSPE(dirname = dirname, countmatfpattern = count_pattern, metadatafpattern = meta_pattern)
spe <- scuttle::addPerCellQC(spe, subsets=list(NegProb = grep("^NegPrb", rownames(spe))))
names(colData(spe)@listData)
## [1] "fov" "cell_ID"
## [3] "Area" "AspectRatio"
## [5] "CenterX_local_px" "CenterY_local_px"
## [7] "Width" "Height"
## [9] "Mean.PanCK" "Max.PanCK"
## [11] "Mean.CD68" "Max.CD68"
## [13] "Mean.MembraneStain_B2M" "Max.MembraneStain_B2M"
## [15] "Mean.CD45" "Max.CD45"
## [17] "Mean.DAPI" "Max.DAPI"
## [19] "sample_id" "sum"
## [21] "detected" "subsets_NegProb_sum"
## [23] "subsets_NegProb_detected" "subsets_NegProb_percent"
## [25] "total"
dim(spe)
## [1] 1010 59284
Selecting only numeric and integer variables from metadata then only variables of interest from the list
temp <- unlist(lapply(spe@colData@listData, function(x) is.numeric(x) | is.integer(x)))
num_variables <- c()
for (i in 1:length(temp)){
num_variables[i] <- ifelse(temp[i] == TRUE, names(temp[i]), NA)
num_variables <- num_variables[!is.na(num_variables)]
}
num_variables
## [1] "fov" "cell_ID"
## [3] "Area" "AspectRatio"
## [5] "CenterX_local_px" "CenterY_local_px"
## [7] "Width" "Height"
## [9] "Mean.PanCK" "Max.PanCK"
## [11] "Mean.CD68" "Max.CD68"
## [13] "Mean.MembraneStain_B2M" "Max.MembraneStain_B2M"
## [15] "Mean.CD45" "Max.CD45"
## [17] "Mean.DAPI" "Max.DAPI"
## [19] "sum" "detected"
## [21] "subsets_NegProb_sum" "subsets_NegProb_detected"
## [23] "subsets_NegProb_percent" "total"
variables <- c("total", "detected")
Approaching dataset preprocessing starting from global-scale, proceeding to local-scale, then fov-focused exploratory analysis
Sample map with fov grid
fov_pattern <- "Run5810_Case2_fov_positions_file.csv"
fovpos_file <- file.path(dirname, list.files(dirname, fov_pattern))
metadata <- data.table::fread(file.path(dirname, list.files(dirname, meta_pattern)))
fov_positions <- data.table::fread(fovpos_file, header = T)
fov_positions <- fov_positions[fov_positions$fov%in%unique(metadata$fov),]
# image theme
my_image_theme <- function(back.color="black",back.border=NA,title.col="white") {
theme(panel.border = element_blank(),
legend.key = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank(),
axis.text.y = element_blank(),
axis.text.x = element_blank(),
panel.grid = element_blank(),
panel.grid.minor = element_blank(),
panel.grid.major = element_blank(),
panel.background = element_rect(fill = "transparent",colour = NA),
plot.title = element_text(color = title.col,hjust = 0.5,face = "bold"),
plot.background = element_rect(fill = "transparent",colour = back.border))
}
fov_size <- 4256
ggplot() +
geom_point(data = metadata, mapping = aes(x = CenterX_global_px, y = CenterY_global_px), colour = "darkmagenta", fill="darkmagenta", size = 0.05, alpha = 0.2) +
annotate("rect", xmin = fov_positions$x_global_px,
xmax = fov_positions$x_global_px + fov_size,
ymin = fov_positions$y_global_px,
ymax = fov_positions$y_global_px + fov_size,
alpha = .2, color = "black",linewidth = 0.2) +
geom_text(aes(x = fov_positions$x_global_px+fov_size/2,
y = fov_positions$y_global_px+fov_size/2,
label = fov_positions$fov), color="black", fontface = "bold") +
ggtitle(sample_name) +
my_image_theme(back.color = "white", back.border = "white", title.col = "black")
Numeric/integer variable dependence on global position
Importing global coordinates from spatialCoords slot since they are not available in colData
colData(spe)$CenterX_global_px <- spatialCoords(spe)[,1]
colData(spe)$CenterY_global_px <- spatialCoords(spe)[,2]
Global coordinate scatterplot colored by variables showing spatial pattern
spe$CenterX_global_um <- spe$CenterX_global_px*0.18
spe$CenterY_global_um <- spe$CenterY_global_px*0.18
spe$logtot <- log10(spe$total)
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "total", point_size = 0.1) + ggtitle("Total counts per cell")
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "logtot", point_size = 0.1) + ggtitle("Log-total counts per cell")
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "detected", point_size = 0.1) + ggtitle("Detected features per cell")
spe$um_area <- spe$Area*0.0324
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "um_area", point_size = 0.1) + ggtitle("Area in μm per cell")
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "AspectRatio", point_size = 0.1) + ggtitle("Width/Height ratio per cell")
spe$target_counts <- spe$total - spe$subsets_NegProb_sum
spe$target_detected <- spe$detected - spe$subsets_NegProb_detected
spe$tot_det_ratio <- spe$target_counts/spe$target_detected
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "tot_det_ratio", point_size = 0.1) + ggtitle("Target counts/detected features ratio per cell")
plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "subsets_NegProb_sum", point_size = 0.1) + ggtitle("Negative control probe counts/detected ratio per cell")
Global coordinate scatterplot colored by binary color code: in black are the cells below the set variable threshold, in grey all the cells above
threshold <- 0.10
spe$var_color <- as.factor(ifelse(colData(spe)$total < quantile(colData(spe)$total, probs = threshold), "black", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Counts below 10% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.05
spe$var_color <- as.factor(ifelse(colData(spe)$total < quantile(colData(spe)$total, probs = threshold), "black", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Counts below 5% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
threshold <- 0.10
spe$var_color <- as.factor(ifelse(colData(spe)$detected < quantile(colData(spe)$detected, probs = threshold), "black", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Detected features below 10% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.05
spe$var_color <- as.factor(ifelse(colData(spe)$detected < quantile(colData(spe)$detected, probs = threshold), "black", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Detected features below 5% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
spe$area_um <- spe$Area*0.0324
threshold <- 0.10
spe$var_color <- as.factor(ifelse(colData(spe)$area_um < quantile(colData(spe)$area_um, probs = threshold), "black", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Area in μm below 10%") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.05
spe$var_color <- as.factor(ifelse(colData(spe)$area_um < quantile(colData(spe)$area_um, probs = threshold), "black", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Area in μm below 5%") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
threshold <- 0.90
spe$var_color <- as.factor(ifelse(colData(spe)$AspectRatio > quantile(colData(spe)$area_um, probs = threshold), "red", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(red = "red",lightgrey = "lightgrey")) + ggtitle("Area in μm above 90% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.95
spe$var_color <- as.factor(ifelse(colData(spe)$AspectRatio > quantile(colData(spe)$area_um, probs = threshold), "red", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(red = "red",lightgrey = "lightgrey")) + ggtitle("Area in μm above 95% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
threshold <- 0.10
spe$var_color <- as.factor(ifelse(colData(spe)$AspectRatio < quantile(colData(spe)$AspectRatio, probs = threshold), "black", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Width/height below 10% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.05
spe$var_color <- as.factor(ifelse(colData(spe)$AspectRatio < quantile(colData(spe)$AspectRatio, probs = threshold), "black", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Width/height below 5% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
threshold <- 0.90
spe$var_color <- as.factor(ifelse(colData(spe)$AspectRatio > quantile(colData(spe)$AspectRatio, probs = threshold), "red", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(red = "red",lightgrey = "lightgrey")) + ggtitle("Width/height above 90% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.95
spe$var_color <- as.factor(ifelse(colData(spe)$AspectRatio > quantile(colData(spe)$AspectRatio, probs = threshold), "red", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(red = "red",lightgrey = "lightgrey")) + ggtitle("Width/height above 95% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
threshold <- 0.10
spe$var_color <- as.factor(ifelse(colData(spe)$tot_det_ratio < quantile(colData(spe)$tot_det_ratio, probs = threshold, na.rm = T), "black", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Target counts/features below 10% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.05
spe$var_color <- as.factor(ifelse(colData(spe)$tot_det_ratio < quantile(colData(spe)$tot_det_ratio, probs = threshold, na.rm = T), "black", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(black = "black",lightgrey = "lightgrey")) + ggtitle("Target counts/features below 5% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
threshold <- 0.50
spe$var_color <- as.factor(ifelse(colData(spe)$subsets_NegProb_sum > quantile(colData(spe)$subsets_NegProb_sum, probs = threshold, na.rm = T), "red", "lightgrey"))
p1 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(red = "red",lightgrey = "lightgrey")) + ggtitle("Negative probe counts above 50% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
threshold <- 0.75
spe$var_color <- as.factor(ifelse(colData(spe)$subsets_NegProb_sum > quantile(colData(spe)$subsets_NegProb_sum, probs = threshold, na.rm = T), "red", "lightgrey"))
p2 <- plotColData(spe, "CenterY_global_um", "CenterX_global_um", colour_by = "var_color",
point_size = 0.1) + scale_color_manual(values = c(red = "red",lightgrey = "lightgrey")) + ggtitle("Negative probe counts above 75% threshold") + theme(legend.position="none")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
plot_grid(p1, p2, ncol = 2)
Heatmap with fovs colored by mean/median values of variable of interest
meta_df <- data.frame(colData(spe))
fov_cellcounts <- as.integer(table(colData(spe)$fov))
meta_df <- group_by(meta_df, fov) |> summarize(mean_total = mean(total), mean_detected = mean(detected), mean_areaum = mean(area_um), mean_aspectratio = mean(AspectRatio), mean_tot_det_ratio = mean(tot_det_ratio), mean_subsets_NegProb_sum = mean(subsets_NegProb_sum), median_total = median(total), median_detected = median(detected), median_areaum = median(area_um),median_aspectratio = median(AspectRatio), median_tot_det_ratio = median(tot_det_ratio), median_subsets_NegProb_sum = median(subsets_NegProb_sum))
meta_df <- cbind(meta_df, fov_cellcounts)
fov_positions <- left_join(fov_positions, meta_df, by = join_by(fov))
fov_size <- 4256
poly_df <- transmute(fov_positions, x_1 = x_global_px, y_1 = y_global_px,
x_2 = x_global_px + fov_size, y_2 = y_global_px,
x_3 = x_global_px + fov_size, y_3 = y_global_px + fov_size,
x_4 = x_global_px, y_4 = y_global_px + fov_size)
str(poly_df)
## Classes 'data.table' and 'data.frame': 45 obs. of 8 variables:
## $ x_1: num 151626 155892 155892 151626 137785 ...
## $ y_1: num 13820 13820 18086 18086 23488 ...
## $ x_2: num 155882 160148 160148 155882 142041 ...
## $ y_2: num 13820 13820 18086 18086 23488 ...
## $ x_3: num 155882 160148 160148 155882 142041 ...
## $ y_3: num 18076 18076 22342 22342 27744 ...
## $ x_4: num 151626 155892 155892 151626 137785 ...
## $ y_4: num 18076 18076 22342 22342 27744 ...
## - attr(*, ".internal.selfref")=<externalptr>
ls_list <- list()
for (k in 1:dim(fov_positions)[1]){
ls_list[[k]] <- st_linestring(matrix(as.numeric(poly_df[k]), ncol = 2, byrow = T))
}
sf_data <- st_as_sf(fov_positions, coords = c("x_global_px", "y_global_px"), geometry = st_sfc(ls_list))
sf_data <- st_cast(sf_data, "POLYGON")
# plot(sf_data)
names(sf_data)
## [1] "fov" "mean_total"
## [3] "mean_detected" "mean_areaum"
## [5] "mean_aspectratio" "mean_tot_det_ratio"
## [7] "mean_subsets_NegProb_sum" "median_total"
## [9] "median_detected" "median_areaum"
## [11] "median_aspectratio" "median_tot_det_ratio"
## [13] "median_subsets_NegProb_sum" "fov_cellcounts"
## [15] "geometry"
sf_data$fov <- as.factor(sf_data$fov)
heatmap_variables <- names(sf_data)[-c(1, length(names(sf_data)), length(names(sf_data))-1)]
for (var in 1:length(heatmap_variables)){
heat_var <- heatmap_variables[var]
p <- ggplot() +
geom_sf(data = sf_data, aes(fill = .data[[heat_var]])) +
scale_fill_viridis_c() +
theme(axis.title.x = element_blank(),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
axis.title.y = element_blank(),
axis.text.y=element_blank(),
axis.ticks.y=element_blank(),
panel.grid = element_blank()) +
geom_text(data = fov_positions, aes(x = fov_positions$x_global_px+fov_size/2,
y = fov_positions$y_global_px+fov_size/2,
label = fov_positions$fov), color="white", size = 2) +
ggtitle(paste0(heat_var))
print(p)
}
Boxplot of each of the numeric/integer variable distribution by fov, useful to see spatial patterns in data
variables <- c(variables, "logtot", "area_um", "AspectRatio", "tot_det_ratio", "subsets_NegProb_sum")
for (var in 1:length(variables)){
n <- grep(paste0("^", variables[var], "$"), colnames(spe@colData))
boxplot(spe@colData@listData[[n]] ~ spe$fov, ylab = colnames(spe@colData)[n], cex = 0.2)
title(paste0(colnames(spe@colData)[n], "~fov"))
}
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 1 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 2 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 8 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 11 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 16 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 19 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 20 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 21 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 22 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 31 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 37 is not drawn
## Warning in bplt(at[i], wid = width[i], stats = z$stats[, i], out =
## z$out[z$group == : Outlier (-Inf) in boxplot 41 is not drawn
Per fov cell counts
fov_df <- data.frame("n_cells" = as.integer(table(spe@colData@listData$fov)), "fov" = as.factor(names(table(spe@colData@listData$fov))))
levels(fov_df$fov) <- order(levels(fov_df$fov), seq(1:length(fov_df$fov)))
ggplot(data = fov_df, aes(x = fov, y = n_cells)) +
geom_col(fill = "#B0C4DE") + theme(legend.position ="none") +
ggtitle("Number of cells per fov")
Total/Area scatterplot
cell_df <- as.data.frame(spe@colData@listData)
ggplot(data = cell_df, aes(x = um_area, y = total)) +
geom_point() +
ggtitle("Counts ~ μm area")
ggplot(data = cell_df, aes(x = um_area, y = AspectRatio)) +
geom_point() +
ggtitle("Counts ~ μm area")
Detected/Area scatterplot
ggplot(data = cell_df, aes(x = um_area, y = detected)) +
geom_point() +
ggtitle("Detected features ~ μm area")
AspectRatio/Area scatterplot
ggplot(data = cell_df, aes(x = um_area, y = AspectRatio)) +
geom_point() +
ggtitle("AspectRatio ~ μm area")
Total/Area scatterplot faceted by fov
cell_df <- as.data.frame(spe@colData@listData)
ggplot(data = cell_df, aes(x = um_area, y = total)) +
geom_point(size = 0.1) +
facet_wrap(~fov) +
ggtitle("Total counts cell ~ μm area by fov")
Detected/Area scatterplot faceted by fov
ggplot(data = cell_df, aes(x = um_area, y = detected)) +
geom_point(size = 0.1) +
facet_wrap(~fov) +
ggtitle("Detected features ~ μm area by fov")
AspectRatio/Area scatterplot faceted by fov
ggplot(data = cell_df, aes(x = um_area, y = AspectRatio)) +
geom_point(size = 0.1) +
facet_wrap(~fov) +
ggtitle("AspectRatio ~ μm area by fov")